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We study the competition between an instantaneous local Coulomb repulsion and a boson medi- 
ated retarded attraction, as described by the Hubbard-Holstein model. Restricting to the case of 
half filling, the ground-state phase diagram and the transitions from antiferromagnetically ordered 
states to charge ordered states are analyzed. The calculations are based on the model in large di- 
mensions, so that dynamical mean field theory can be applied, and the associated impurity problem 
is solved using the numerical renormalization group method. The transition is found to occur when 
electron-electron coupling strength U and the induced interaction A due to electron-phonon coupling 
approximately coincide, U ~ A. We find a continuous transition for small coupling and large u>o, 
and a discontinuous one for large coupling and/or small ojo- We present results for the order param- 
eters, the static expectation values for the electrons and phonons, and the corresponding spectral 
functions. They illustrate the different types of behavior to be seen near the transitions. Addi- 
tionally, the quasiparticle properties are calculated in the normal state, which leads to a consistent 
interpretation of the low energy excitations. 

PACS numbers: 71.10.Fd, 71.27.+a,71.30.+h,75.20.-g, 71. 10. Ay 



I. INTRODUCTION 

A feature of strongly correlated systems is the existence 
of competing interactions on low energy scales which can 
lead to different types of symmetry breaking and differ- 
ent ground states. There can be various forms of mag- 
netic order, superconducting or charge ordered states; 
there may also be transitions between these states, and 
in some cases they even coexist. For instance, compounds 
such as the vanaditesi, high T c cuprates 2 , fullerides 3 , 
manganites^ and organic salts^ possess rather involved 
phase diagrams and, to understand them, an analysis of 
the competition between the different interactions will be 
important. 

Here we study the competing effects between an in- 
stantaneous local Coulomb repulsion and the retarded 
interaction induced by a coupling to an optical phonon 
mode using the Hubbard-Holstein (HH) model. We con- 
sider the competition between two types of order, anti- 
ferromagnetic (AFM) and charge order (CO), which can 
occur in the model at half filling. The emphasis will be on 
treating the phonons fully quantum mechanically and in 
allowing for arbitrary coupling strengths, so that the full 
interaction parameter regime can be investigated. This 
is possible if we use the infinite dimensional version of 
Hubbard-Holstein model so that we can apply the dy- 
namical mean field theory (DMFT), which becomes ex- 
act in this limit. The numerical renormalization group 
(NRG) method is then used to solve the associated effec- 
tive impurity problem. This permits one to handle both 
strong electron-electron and strong electron-phonon in- 
teractions as well as a wide range of phonon frequencies. 
We focus on the ground state and spectral properties of 
electrons and phonons at zero temperature. 

The infinite dimensional Hubbard. Holstein and com- 



bined HH models have received considerable attention 
in the pas&~—. For the pure Holstein case, Freericks 
et al£& found instabilities to charge order and super- 
conductivity by quantum Monte Carlo (QMC) and it- 
erated perturbation theory for different filling factors. 
At half filling, Benedetti and Zeyher— , and Hague 
and D'Ambrumenil 2 ^, investigated the normal state and 
found a breakdown of Migdal-Eliashberg theory when a 
lattice instability develops for stronger electron-phonon 
coupling. The charge ordered ground state and phase di- 
agram in the adiabatic limit has been analyzed by Ciuchi 
et alfii. It was shown there that the weak and strong 
coupling CO states are smoothly connected. 

For the Hubbard-Holstein model in the absence of 
long range order, the phase diagram of the paramagnetic 
(PM), bipolaronic (BP) phases and the metal-insulator 
(MI) transition has been establishe d 16 : 24 . Another re- 
cent study of the model without long range order deals 
with the topic of polaron formatio n 18 : 19 with finite elec- 
tron density, extending the original work of Holstein 25 
who considered the single electron case only. The occur- 
rence of superconductivity was studied in Ref. 

mum 

In a two site calculation Takada found superconductivity 
with off-site pairing at half filling in a very small param- 
eter regime in the antiadiabatic region 21 . The effect of 
phonons on the quasiparticle excitations in the presence 
of AFM has also been investigated 26 . There have also 
been extensive calculations for the one dimensional ver- 
sion of model2i" ~— , which we shall comment on briefly 
later. 

Our analysis here of the HH model will extend the ear- 
lier work by allowing for AFM and CO states, which are 
the dominant instabilities at half filling. We study the 
transitions between these states. This will give a more 
complete picture of the phase diagram and the proper- 



2 



ties of the model in the ordered phases. In the regions of 
the phase diagram with CO, we also obtained supercon- 
ducting solutions, but the CO states were found to have 
lower energy. We calculate the static and dynamic prop- 
erties in both these types of broken symmetry phases. 
The paper is structured as follows. In Sec. II we specify 
the formal setup of the HH model and the DMFT-NRG 
method. We also give explicit expressions for the differ- 
ent contributions to the total energy. The dependence 
of these on the interactions is discussed in detail in Sec. 
IV. Before that in Sec. Ill, we discuss the global phase 
diagram, the order parameters and static quantities and 
their dependence on U, A and ljq. Sec. V explores the 
normal state properties of the HH model which helps to 
understand the ground state phase diagram and transi- 
tion. In Sec. VI we discuss how the bosonic properties 
are modified by the coupling to the electronic system. In 
Sec. VII we present results for the electronic and bosonic 
spectral functions, before concluding in Sec. VIII. 



II. MODEL AND DMFT-NRG SETUP 



The Hamiltonian for the HH model is given by 



(1) 



+UJQ 



E^+sE^+^KE^- 1 



cj a creates an electron at lattice site % with spin a, and 

b[ a phonon with oscillator frequency uiq, ni_ a — c\ a Ci iCT . 
The electrons interact locally with strength U, and their 
density is coupled to an optical phonon mode with cou- 
pling constant g. We have set the ionic mass to M = 1 
in (Q]). The local oscillator displacement is related to the 
bosonic operators by X{ = (6, + bj)/ \/2ojq. where h = 1, 
and one can define a characteristic length xq = 1/ ^/ujq 
for the oscillator. In Appendix lAl we give the details for 
a mean field calculation in the adiabatic limit for this 
model. 

For our calculations we assume a bipartite lattice with 
A and B sublattice, where the matrix Green's function 
can be written in the form 



1 



'CA,cr(w)Cfl,a(w)-efcV £k ^A^\ 



(2) 



with i^ airT (uj) = ui + \x a ,o — S Q:(T (aj), a = A,B, and 
fc-independent self-energy 3? . For commensurate charge 
order we have pA.a = P — h c , ps.a = p> + h c and 
£b j0 -(w) = Un - Ex,(r(-w)*, with n = (n A + n B )/2, 
n a = ^2 a n a,tn n u,cr = {n a ,a) ■ For the AFM order one 
has fiA.a = f 1 ~ &h s , PB,a = p + <jh s , and the condition 
Eb ](T (uj) = Si We consider solutions of exclusive 

AFM or CO, where the symmetry breaking fields vanish, 
h c , h s — > 0. 



In the case with symmetry breaking, the effective Weiss 
field is a 2 x 2 matrix Q_~ 1 (t) . The DMFT self-consistency 
equation in this case reads^ 3 . 



(3) 



The matrix of local lattice Green's functions G CT (w) — 
V^Sfc Q-k ffC^O i s obtained by integrating over the den- 
sity of states, l/iVY^fc /(fffe) = Jde po(e)f(e). We as- 



sume a semi-elliptic DOS, po(e) = 2\/ D 2 — e 2 /ttD 2 cor- 
responding to a Bethe lattice in all the following calcula- 
tions. In the DMFT this local Green's function, and the 
self-energy are identified with the corresponding quan- 
tities for an effective impurity model^ 3 .. One focuses for 
the calculations on the properties of the A-sublattice. We 
can take the form of this impurity model to correspond to 
an Anderson-Holstein impurity model^i and calculations 
are carried out as detailed, for instance in Ref. 35.3(| We 
solve the effective impurity problem with the numerical 
renormalization group 37 ! 38 (NRG) adapted to these cases 
with symmetry breaking. The NRG has been shown to 
be very successful for calculating the local dynamic re- 
sponse functions, and we use the recent approac h 39 ! 40 
based on the complete basis set proposed by Anders and 
Schiller—. For the logarithmic discretization parameter 
we take the value A = 1.8 and keep about 1000 states at 
each iteration. The bosonic Hilbert space is restricted to 
a maximum of 50 states. 

In the AFM case the A-sublattice magnetization, 
'I'afm = m A = ( n A,-\ — nA.i)/2 serves as an order pa- 
rameter. For CO we define $ co = (tia — l)/2. 

To find the ground state of the system we calculate 
the ground state energy E to t = (H) /N of the HH Hamil- 
tonian (TTJ) in the different phases. This gives generally, 



Etot = -Ekin + En + E g + E p h- 



^ T- -pn- (4) 

The first term is the kinetic energy, which reads 



-Ekin = E j dSk /°o( £ fc) £ fc J^U) f(u)pAB,k,<r{u), (5) 

where pAB,k,a(u) = -lmGAB,h,a{^)/^ for the off- 
diagonal Green's function in $2$ and f(u>) is the Fermi 
function, where f(uj) = 9{—lj) at zero temperature. In 
the non-interacting case this expression can be evalu- 
ated analytically and we find for half filling, p = 0, 
£° n = -4L>/3tt, which for D = 2 is E° in ~ -0.849. 
This can be used as reference energy. More specifically 
one finds 

E k in = E J d£k Po( £ k)£k Jduj f((jj)g k ,a(uj), (6) 



where 



9k ,<r(ui) = Im 



(7) 
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The interaction energies Ejj . E g can be calculated from 
expectation values, 

E v = \ $><*.t<+>' E g = | ^<(6 a + &£)(n a - 1)> 

a a 

We distinguish between A- and _B-sublattice values, 
which are equal in the AFM case, but not for the CO 
case. There we use the operator identity fiB,a — 1 — fiA,a 
(particle hole transformation), at half filling, which yields 
{fiB^fiB^) = («A,t?W) + 1 - («a)- The terms for E g 
turn out to be equal on A and B sublattice as both terms 
will contribute with opposite value. 

III. PHASE DIAGRAM AND STATIC 
PROPERTIES OF THE ELECTRONS 

As the electron-phonon coupling in (Q} is linear, the 
bosonic field can be integrated out in a path integral 
framework, which yields a purely electronic theory with 
an effective electron-electron interaction of the form 

U cS (iu) = U + g 2 D (oj), (8) 

where the free phonon propagator on the real axis is 
D°(ui) = 2u>o/(uj 2 — u! 2 ). These are the competing in- 
teractions on different energy scales. There is a limiting 
case wo — > oo, the antiadiabatic limit, where A = 2g 2 /ujq 
is kept fixed. In this limit, U e g(oj) becomes independent 
of uj and tends to U e g = U — A, so that the model then 
becomes equivalent to a Hubbard model with U = U e s. 
Generally, the situation is more complicated. For large 
uj the Coulomb repulsion U is dominant in ([S} as D°(ui) 
goes to zero. However, ujq enters as a relevant energy 
scale at lower energy and for \uj\ < ujq the competition 
between the bare interactions is most important. The 
DMFT calculations deal with these competing interac- 
tions on different energy scales, and as a result the ground 
state phase diagram of the infinite dimensional HH model 
at half filling emerges. This phase diagram was already 
presented earlier—. In order to give a comprehensive ac- 
count in this paper, we include some of these results in 
the following. 

We comment briefly on the notation and energy 
scales. In electron-phonon physics, specifically in Migdal- 
Eliashberg theory, the electron phonon coupling is usu- 
ally specified by a dimensionless parameter, often called 
A, involving the coupling strength and an electronic scale. 
In contrast, in this paper A has the dimension of an en- 
ergy, and it is compared with the Coulomb repulsion U. 
The dimensionless parameter for electron-phonon cou- 
pling corresponding to the usual convention would be 
A = po(0)A, where po(0) = 2/(Dw). For physical op- 
tical phonons the scale luq is expected to be roughly 
ujo ~ 0.01I-F — 0.2W, where W is the electronic band- 
width. The purpose of the paper is to characterize 
the electron-phonon system quite generally, with a tun- 
able parameter ujq, which for many results we chose as 



ujq = 0.15W neither too close to the adiabatic nor to the 
antiadiabatic limiting cases for illustration of intermedi- 
ate behavior. It should be noted that the value is large 
for most realistic electron-phonon systems with the pos- 
sible exception of the alkali doped fullerides^ 3 -. For the 
former cases smaller values such as luq = 0.05W could 
serve as a better guideline. 



A. Phase diagram 

We first present the complete phase diagram as shown 
in Fig. [TJ The overall energy scale is set by the band- 
width W — 4t — 4, and the phonon frequency ujq = 0.6t 
was chosen. Many results such as the global phase dia- 
gram arc similar for different choices of ujq, but we will 
also point out the differences appearing for other values 
of ujq- The limiting cases of the phase diagram can be un- 
derstood on a qualitative level. The U-axis corresponds 
to the pure repulsive Hubbard model which is known to 
be AFM ordered at weak coupling^ for a bipartite lat- 
tice, and this order is smoothly connected to the strong 
coupling Heisenberg AF M 35 i 36 . Also along the A-axis, the 
corresponding pure Holstein model has a charge ordered 
ground state for g > and ujq > 0, such that the limits of 
weak and strong coupling are smoothly connected^. For 
finite U and g we find that the transition line is close to 
the line A = U with a small tendency towards A > U. 
For larger values of the interactions we have included a 
dashed line (o) above which our DMFT-NRG calcula- 
tions find solutions with finite $ co and another one (o) 
below which ^afm is finite. 

For weaker coupling (U < 3) the order parameters be- 
come very small close to the line U B s = 0. There are 
strong indications that the transition proceeds both di- 
rectly and continuously from one type of ordered to the 
other ordered state, such that only at the critical point 
both order parameters vanish. For instance, we find that 
all the relevant response quantities behave continuously. 
A direct order to order transition is found in the antia- 
diabatic limit, where for any finite J7 G ff > the system is 
AFM ordered^ and for U cS < in the CO or SC stated. 
We have therefore a continuous transition from an or- 
dered to an ordered state with vanishing order parame- 
ters at the transition. There is no reason to expect that 
the transition becomes discontinuous immediately when 
wo is decreased from infinity to a finite value ujq 3> W. 
In addition, in the DMFT-NRG calculations we find that 
near the transitio, the smaller o>o is, the larger the order 
parametern becomes, as detailed later in Fig. [7] Thus 
we conclude that the direct continuous transition sce- 
nario persists for weak coupling and finite ujq. Mean field 
calculations (see Appendix [XJ in the adiabatic limit sup- 
port the picture of a direct transition from one ordered 
state to the other, however, the transition is always dis- 
continuous then. Also the local effective quasiparticle 
interaction U r (which will be discussed more fully later) 
is observed to change sign at J7 G ff = 0, which is consistent 
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We see that only for large coupling do the phase bound- 
aries merge, while for smaller couplings other scales are 
important. For the case ojo = 0.2 this has been ana- 
lyzed earlie r 16 i 17 ' 24 . The Mott transition, as obtained on 
increasing U for fixed g is only little affected by the addi- 
tional electron-phonon coupling, which is manifested in 
a shift of the critical U c for the transition. The metal- 
bipolaron transition, observed when increasing g for fixed 
U, is of second order for smaller interactions and becomes 
of first order for large interactions, and thus similarities 
with the ground state behavior are found. Note that 
the CO state and the bipolaronic state are different, as 
in the latter no symmetry is broken and the occupation 
expectation value is always 1, whereas in the CO state 
n A ^ n B . For details we refer to Refs. [l6l fl7l|24l 



FIG. 1: (Color online) The phase diagram of antiferromag- 
netic (AFM) and charge order (CO) in the U- A- plane. The 
thin black line A ~ U gives a continuous transition and the 
thick line a discontinuous one. A non-zero order parameter 
was found above dashed line with diamonds for CO and below 
the dashed line with circles for AFM order. The dashed lines 
with points gives the transition for phases with no long range 
order, a paramagnetic metallic (PM), bipolaronic (BP) and 
Mott insulator (MI). 



with a change of ground state there. 

For larger couplings we have a parameter regime where 
we can find finite $ a f m and $ co . In this regime the tran- 
sition from the AFM to the CO turns out to be discontin- 
uous and one can identify a point on the transition where 
the nature of transition changes. The calculation of the 
total ground state energy (see inset of Fig. HJ) shows that 
also here the transition occurs approximately at U c s ~ 0. 
The behavior of the model along the line U = A has been 
studied in Ref4£. 

The HH model in one dimension has been studied in 
great detail with efficient numerical methods 27-29 . One 
finds a Mott insulator with strong antiferromagnetic cor- 
relations, but no long range order, when [/off > 0, and a 
Peierls charge density wave (CDW) insulator for U c s < 0. 
There is, however, a metallic region with finite spin gap, 
but no charge gap in between these two phases. The 
transition line to the CDW state appears for values of 
A a bit larger than U, similar to what is observed in 
our calculations but more pronounced. For larger U this 
intermediate region shrinks until we get a direct first or- 
der Mott-Peierls transition. A major difference with the 
high dimensional results is the real symmetry breaking 
in our case as well as the existence of the intermediate 
region, for which we find no indication here. Neverthe- 
less the mentioned similarities and also recent results in 
d = 2 (adiabatic limit) 47 , suggest that the general features 
of the phase diagram in Fig. [T] might be quite general, 
largely independent of dimensionality. 

In Fig. [T]we have also included the phase boundaries 
of the HH model when no long range order is allowed for. 



B. Order parameter as a function of A 

We consider here the way the two types of order pa- 
rameter change as a function of A for fixed values of U 
both in the weak coupling and strong coupling regime. 
A weaker coupling case with U — 2 is shown in Fig. 
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FIG. 2: (Color online) The expectation values <3? for U — 2 as 
a function of A. The inset shows the total energy. 

We can see that the AFM order is decreased when the 
electron phonon-coupling is increased, as the repulsion is 
reduced. Near A = U the ordering scale is very small 
(< I(T 3 ) and cannot be resolved in our DMFT-NRG 
calculations. For A > U the $ co shows a steep rise with 
A. For this weak coupling case we can study the limit 
ljq — >• and compare with the corresponding static mean 
field theory (details see appendix [A"|l . For U = 2 the 
numerical results are shown in Fig. [3J 
The solutions are for situations where the order is exclu- 
sive, i.e. only one of the order parameters is nonzero. In 
the mean field calculation the AFM order parameter is 
larger than in the DMFT case, and when A is increased 
is seen not to be affected by the electron phonon cou- 
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FIG. 3: (Color online) Mean field result for the order param- 
eter and total energy shown in the inset. 



FIG. 4: (Color online) The expectation values $ for U = 5 as 
a function of A. 



pling as long as $ co = 0. In contrast the CO order state 
<E> co feels the {/-term, but increases with A. For U = A 
the mean field equations give order parameters that coin- 
cide, and once 4> co exceeds $ a fm the charge ordered state 
possesses the lowest energy as can be seen from equation 
(IA13|) and the inset in Fig. [3] . We can infer from this 
weak coupling result that there is a direct, discontinuous 
transition from an ordered to an ordered state at U = A 
for uj — > 0. As seen in Fig. [5] the behavior is strongly 
modified when quantum fluctuations are included and ujq 
is well finite, as the order parameters are influenced much 
more by the presence of the competing interaction. It is 
numerically not possible to study the limit ujq — > within 
our DMFT approach due to the increase in the bosonic 
Hilbert space. The analysis for different values of ujq is 
however conform with the trends discussed here (see Fig. 

ED. 

For larger couplings we saw that the two dashed lines 
in Fig. [1] cross, which means that we have a parameter 
regime where we find finite $ a f m and $ co . An example 
for this behavior is shown in Fig. @]for U = 5. 
The transition here is seen to be rather sharp. The calcu- 
lation of the total ground state energy (see inset) shows 
that the transition occurs approximately at U c g = 0, in 
fact it occurs for small negative U e g, i.e., on the A > U 
side. A number of quantities such as the double occu- 
pancy (n^hi) (see section UlI CI) show discontinuities at 
the transition. The total energy (see inset of Fig. [3} is 
continuous function of A. It shows, however, a kink at 
the transition, such that a first derivatives will be discon- 
tinuous. 



C. Double occupancy 

A characteristic quantity for the electronic part of the 
system is the expectation value for the local double occu- 
pancy. This is a homogeneous quantity in the normal (N) 



and AFM phase. In the charge ordered phase it differs for 
A- and £?-sublattice, and (hA^fiA,i) and {n 3,^3,1.) are 
given as detailed above when discussing the energy due 
to the Hubbard interaction term. We can compare the 
quantities by taking the average over the two sublattices, 
(n-|-nj,) = E a ("a,t"a,|}/^ I n Fig- El we show the results 
for U = 2, 5 as a function of A. We have included the N, 
AFM and CO state, and for the latter, the averaged as 
well as the sublattice quantities. 

For the continuous transition at U = 2 we find a small 
nearly linear increase with A for the N and AFM state. 
The values for N and AFM state are very similar with 
the double occupancy being larger in the normal phase, 
where the charge carriers are more mobile than in the 
ordered phase. This is different from the strong coupling 
case, where the N state is in the Mott phase, and the 
double occupancy is lower than in the AFM case. This 
is due to the super-exchange mechanism which leads to 
a kinetic energy gain and a slightly higher mobility and 
also double occupancy. When A > U and charge order 
sets in, the double occupancy shows a steep rise on the 
A-sublattice and decrease on the i?-sublattice. As can be 
seen in the top panel of Fig. [5l the A-dependence is con- 
tinuous for U = 2. At approximately A = 3 the N state 
metal-bipolaron transition occurs, where (fifhi) in the 
normal state increases rapidly (but also continuously). 
(n-|-nj,) is then larger than in the CO state which can be 
understood by thinking about the fact that the CO state 
wins energetically against the BP state through kinetic 
energy due to pair hopping (see Sec. IIVI) . 

The overall behavior in the lower panel of Fig. [5] for 
U = 5 is similar, (n-j-n^ increases slightly with A for a cer- 
tain range and more rapidly near the transition. There, 
we can clearly see the discontinuous behavior for A ~ 5 
and (n^ni) jumps from the value 0.09 in the AFM state 
to 0.42 in the CO state. The BP state again has a larger 
double occupancy than the average of the CO state. 
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FIG. 5: (Color online) The expectation values (n^ni } for 
(7 = 2 (upper panel) and U — 5 (lower panel) as a function 
of A for the N, AFM and CO state. 



FIG. 6: (Color online) The expectation values $ for A ~ 2.13 
(top) and $ for A ~ 4.8 (bottom) as a function of U. The 
insets show the total ground state energy. 



D. Order parameter as a function of U 

So far we have analyzed the transition as a function of 
A. Similarly, we can look at the properties of the model 
as a function of the on-site repulsion U. We focus on the 
order parameters and look at the cases of a fixed A ~ 2.13 
(g = 0.8) and A ~ 4.8 (g = 1.2) shown in Fig. 1 
In both cases for small U the CO state dominates, but 
on increasing U, $ co is driven to zero and $ a fm becomes 
finite. Similar to the weaker coupling case for fixed U for 
A ~ 2.13, the order parameters approach zero at the tran- 
sition. We also have continuous transition here, which is 
visible in the total energy plotted as an inset. 

In contrast for larger A ~ 4.8 (lower panel in Fig. [6]) 
the picture is as for the larger coupling cases above, where 
the transition occurs with finite order parameters, which 
change discontinuously. Also here, as seen in the inset, 
the total ground state energy displays a kink at the tran- 
sition. Summarizing, we can say that the transition oc- 
curs in similar fashion as a function of U or A. Depending 
on the magnitude of the coupling constants continuous or 



discontinuous transitions can be observed. 



E. Dependence of the order parameter on uiq 

So far we have found that the transition from the CO 
to the AFM state occurs close to U e s = independently 
of the other parameters. For fixed ujq, when studying 
the order parameters or double occupancy, one could see 
quite a different dependence on A or U near the transition 
with a continuous and discontinuous behavior. In the fol- 
lowing we study how the order parameters depend on U e g 
near the transition for different values of wq. We illus- 
trate this in Fig. [71 where we plot the respective order 
parameters for different ujq. including the case ujq — > oo, 
which is given by the pure Hubbard model with local in- 
teraction U c s- We have held U — 3 fixed and varied g to 
obtain the desired values. 

As one can see clearly the "sharpness" of the transition 
increases when ujq is decreased. For ujq — > oo the order 
parameters approach zero in a similar exponential form 
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FIG. 7: (Color online) The expectation values $ for a range 
of effective U e s = — 2, 2 for fixed [7 = 3 varying g. We show 
the results for ojq = 0.2 (full line), uuq — 0.6 (dashed line), and 
cjo — > oo (dot-dashed line). 

as in mean field theory^. The form of the transition is 
then symmetric with respect to U e g — 0. For finite loq 
when tuning £7 e g with g the transition will be asymmetric 
and the order parameter is decreased much less when the 
effective interaction is close to 0, as seen most pronounced 
for the case luq — 0.2 (full line). In the adiabatic limit, 
cjo — > 0, we always expect a discontinuous transition at 
zero temperature. For finite luq we have the competition 
of the Hubbard repulsion and phonon-induced attraction. 
The effect of latter enters at lower energies for smaller 
values of luq. This might explain why the AFM order is 
more stable then. To establish CO as the ground state, it 
seems to be mainly necessary that £/ < A, and the large U 
at higher energies does not spoil this. Retardation effects 
seem to play hardly any role at half filling for these static 
orders. 

We can conclude from this section that with general- 
ity the AFM-CO quantum phase transition occurs ap- 
proximately when electron interaction parameter U and 
phonon attraction A are equal. The behavior near the 
transition, e.g. the order parameter depends, however, 
very much on the the interaction strength as well as the 
phonon frequencies. Small phonon frequencies and large 
interactions lead to discontinuous behavior, whereas for 
large phonon frequencies the competing interactions lead 
to more cancellations, reduced order and much evidence 
for continuous transitions. Then there exists a point on 
the transition line At c — U, separating continuous and 
discontinuous transitions. We find that the value At c in- 
creases with ujq. In the limiting case of luq — > oo, there are 
only continuous transitions such that At c = oo, whereas 
we expect that for ujq — > there are only discontinu- 
ous transitions and At c = 0. It would be of interest to 
explore how At c varies for finite loq as a function of tem- 
perature. As temperature tends to increase fluctuations 
and decrease order, the naive expectation would be that 



Atc increases with temperature. 



IV. DIFFERENT CONTRIBUTIONS TO THE 
TOTAL ENERGY 

We have used the total ground state energy E tot to de- 
cide whether the AFM or CO state is the ground state. 
In this section we give details of the different contribu- 
tions to Etot in (HI) , and their dependence on the electron- 
phonon coupling for fixed U. Let us first remark generally 
on the energy of the ordered state in comparison with the 
normal (N) state. In the half filled pure Hubbard model 
at weak coupling the AFM state has lower potential en- 
ergy than the N state, but higher kinetic energy. At 
strong coupling the AFM has lower kinetic energy than 
the N state (exchange term) , but higher potential energy 
than the N state (Mott insulator). In the pure Holstein 
model the CO state has a lower potential energy at weak 
and intermediate coupling and higher kinetic energy due 
to localization. At strong coupling when the N state is 
insulating (bipolaronic) the energy is lowered in the CO 
state due to lower kinetic energy (pair hopping). 

In the situation with finite U and <?, we have a compe- 
tition between the different terms. The AFM state will 
usually have smaller Ejj, since double occupancy is lower, 
whereas the CO state possesses larger Ejj but contribu- 
tions from Eg lower the energy. However, we also have to 
take into account the contribution 7? p h, which is larger 
in the CO state. 

In Fig. [8] we see the behavior of kinetic energy for the 
electrons E^ n and oscillator energy of the phonons Eph, 
where the energy u>o/2 for zero point motion was omitted. 




FIG. 8: (Color online) Comparison of the kinetic energy 7?ki n 
and phonon energy E p h for U — 2 (left) and (7 = 5 (right) as 
a function of A in CO, AFM and N state. 

We show the quantities for the values U — 2 (left) and 
£7 = 5 (right) as a function of A in CO, AFM and N 



state. As a reference energy recall that for the free sys- 
tem E® in ~ —0.849. Similar to the pure Hubbard or 
Holstein model the electronic kinetic energy of the AFM 
and CO state is larger than the N (metallic) state, whilst 
the energy gain comes from the interaction energies, Eu 
for the AFM state and E g for the CO state, see Fig. 
both for U = 2 and [7 = 5. The situation is different 
for the BP state (A > 3.1 for U = 2 and A > 5.4 for 
[7 = 5). Then the CO state wins through lower kinetic 
energy, whilst the interaction energy is lower in the BP 
state. 

The phonon energy does not change much with A for 
values A < U in the N and AFM state, but increases 
rapidly for A > U both in the N and CO state. This is 
mainly due to the increase in potential energy ujq(x 2 }/2, 
as will be seen in detail in Section |VI] when we discuss 
(x 2 ). 
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to — A. For the CO state this can be understood through 
a factorization of the expectation value in E g in the elec- 
tronic and phonon part. With equation © discussed 
later one finds then E g — > — 4A$ 2 Q . Since in the ordered 
state $ co = 1/2 the behavior follows. 



V. QUASIPARTICLE PROPERTIES IN THE 
NORMAL STATE 

In a large part of the phase diagram in Fig. [TJ near 
the AFM-CO transition, the system without symmetry 
breaking is in the N metallic state, which is a Fermi 
liquid. We can gain insight into the properties of the 
system, when we analyze the quasiparticle properties of 
the N state. The states with broken symmetry can then 
be viewed as instabilities of the Fermi liquid state. The 
renormalization factor z, which is related to the weight 
of the quasiparticle peak and for a fc-indcpcndent self- 
energy to the inverse of the effective mass of the quasi- 
particles, to* /mo — z^ 1 . It can be calculated from the 
derivative of the self-energy as well as from the analysis of 
the NRG low energy excitations at the fixed point. From 
the latter procedure, one can also deduce a local effective 
quasiparticle interaction U r by comparing the energy of 
the lowest two-particle excitation energy E pp with the en- 
ergy of two one-particle excitations E p , U r ~ E pp — 2E p . 
For details we refer to Ref. [3^1491 



A. Renormalized parameters as function of A 

First we discuss how the quasiparticle weight z varies 
with A, when U is held fixed. It is shown for various 
values of U as a function of A in Fig. [TU1 



FIG. 9: (Color online) Comparison of the interaction energies 
Eu and E g for U — 2 (left) and U = 5 (right) as a function 
of A in CO, AFM and N state. 

The energy of the AFM state depends relatively lit- 
tle on the electron-phonon coupling, i.e. if the order is 
not destroyed, the phonons have a minor effect on the 
static properties. Eu has a small value for weak electron- 
phonon coupling and tends to U/2 in the case of strong 
coupling for the BP and CO state. For weaker repulsion, 
U = 2, the energies vary continuously with A. However, 
for the stronger coupling case U = 5 we can clearly see 
the discontinuities at the AFM-CO transition, U ~ A. 
E p h jumps from a low value in the AFM state to a large 
value in the CO state. Due to the increase in double 
occupancy also Ejj increases suddenly at the transitions. 
Both of these energetically unfavorable contributions for 
the CO state are counterbalanced by the abrupt decrease 
of Eg, as the electron-phonon gives strong binding energy 
in the CO state. As a result the total energy is continuous 
at the transition, but shows kink. 

For large A the interaction energy E g is proportional 
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FIG. 10: (Color online) The 2-factor for various U as a func- 
tion of A. 

We can identify different effects of the electron-phonon 
coupling. For U — and small values of U, increasing A 
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leads to polaron formation and localization of the charge 
carriers, which results in a reduction of z and a larger 
effective mass. Eventually, the metal-bipolaron transi- 
tion is reached, where z — > 0. For larger values of U 
the electrons are already renormalized for A = due to 
the Coulomb interaction. The first effect of increasing 
the electron-phonon coupling is to reduce this effect and 
z increases with A. Note that the effect is substantially 
less than what would be expected for a pure Hubbard 
model with U e ff . The maximal value obtained occurs ap- 
proximately when A ~ U, i.e. the renormalization effects 
cancel there to the largest extent leading to a least en- 
hanced effective mass. Near the M-BP transition to* /mo 
diverges as for the weak U case. Apart from the approx- 
imate maximum no particular characteristic behavior is 
seen in z near the AFM-CO transition, U ~ A. 

As shown in Fig. [TT] the effective quasiparticle inter- 
action U r varies between positive and negative values 
depending on U and A. In the cases for finite U it starts 
repulsive and goes to zero approximately where U e g does. 
Then it becomes negative such that there is an effective 
attraction between quasiparticle excitations. There is a 
slight shift towards A > U for this sign change to occur 
for larger values of U. This is in line with the earlier 
observation that the CO state becomes the ground state 
when A > U. 

The change of sign of U r when U ~ A can be related 
to the maximum found for z at this point in Fig 1101 At 
U ~ A, U r ~ 0, the quasiparticles are effectively non- 
interacting. As A is varied from this point, both \U e $\ and 
\U r \ increase which causes a further renormalization of 
the quasiparticles and a reduction in the value of z. The 
decrease in z from U r = occurs irrespective of whether 
U r < 0, as in the pure Holstein model, or U r > 0, as in 
the Hubbard model. 




FIG. 11: The effective quasiparticle interaction U r for various 
U as a function of A. 



Taking the viewpoint of instabilities of the Fermi liq- 
uid, one can infer that the ground state, AFM or CO, 
is determined by whether the low energy quasiparticles 



interact attractively or repulsively. The sign change in 
turn occurs when the bare parameters U , and A on high 
energies are equal. One might have expected that in re- 
ducing the energy scale u down to the phonon frequency 
loq that the main retardation effects would renormalize 
U to some effective value U, where U would be of the or- 
der of U r in the pure Hubbard model, which is such that 
U r <C U. The attractive term induced by the phonons 
would then contribute for to <C ojq ■ The change of sign of 
the quasiparticle interaction would then be expected to 
occur when U ~ A, which would correspond to a much 
smaller value of A than A ~ U. The fact that the transi- 
tion and the change of sign of U r are found to occur when 
U ~ A indicates that U term and the A term are renor- 
malized in a similar way as the energy scale is reduced. 
As noted before, however, when discussing the order pa- 
rameter, the system cannot be described simply by an 
effective Hubbard model and both the electronic interac- 
tion U and A play a role in determining the properties 
of the system in a certain phase. For instance, the CO 
order parameter in Fig. 0] corresponding to U e s ~ —0.01 
is $ co — 0.4, but the result would be close to zero for the 
pure Hubbard model with this interaction on all energy 
scales. 

When A exceeds a certain value \U r \ decreases again, 
and we can also see that not only z but also the effective 
quasiparticle interaction U r goes to zero at the M-BP 
transition. It is of interest to study the combined quan- 
tity U r /z = U r m*, which is plotted in Fig. HU 




FIG. 12: The quantity U r / z for various U as a function of A. 

This product of the effective quasiparticle interaction and 
effective mass takes into account the aspect of the lo- 
calization tendency of the quasiparticles as well as their 
residual interaction. We can see that it shows a more 
universal behavior for the different cases. It decreases 
monotonically as a function of A for the given values of 
U. Close to the M-BP transition in all cases U r m* tends 
to a value between 2.5 and 3. One is therefore tempted to 
identify this as the relevant quantity for an instability of 
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the metallic state, such that the M-BP transition occurs. 



B. Renormalized parameters as function of U 

Similarly, we can also study the behavior of the quasi- 
particle properties for fixed A as a function of U. The 
quasiparticle weight z is shown for various values of A as 
a function of U in Fig. [T51 




FIG. 13: The z-factor for various A as a function of U. 

For A = 0, z is monotonically driven to zero when in- 
creasing U to the Mott transition. For finite values of A 
increasing U leads to an increase in z when we come from 
the bipolaronic state and thus to a de-renormalization as 
\U e g\ decreases. Similar to the case observed above but 
more pronounced, one finds z to be approximately max- 
imal near the AFM-CO transition, U ~ A, i.e. when 
U c s a s well as U r are close to zero. Different renor- 
malization effects have then canceled to maximal ex- 
tent, leaving relatively weakly renormalized, nearly non- 
interacting quasiparticles. Again it is of interest to study 
the quantity U r m*, which is plotted in Fig. [TJ] 
As noted before it goes to zero where C/ e ff does. Then it 
changes sign and increases with U as long as the system is 
in the metallic state. Close to the Mott transition it ap- 
proaches a value 2.5-2.7 which is approximately the same 
for different values of A. It seems therefore that for this 
metal insulator transition a universal value determines 
when it occurs. 



VI. PROPERTIES OF THE PHONONS 

So far we have studied the CO-AFM transition via 
the electronic properties, which where greatly influenced 
by the interaction with the bosonic modes. In turn the 
properties of the local harmonic oscillator modes are also 
modified by the coupling to the electronic degrees of the 
freedom, depending on the different coupling strengths. 
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FIG. 14: The quantity U r / z for various A as a function of U. 



First we consider the phonon number expectation value 
"ph = (b'b). It is expected to increase in the CO (BP) 
state as a high probability of local double occupation 
leads to a charge redistribution and thus a (temporal) dis- 
placement of the lattice ions. This means that phonons 
become excited, multiplied with loq this gives the ener- 
getic contribution E p h = wo^ph This was discussed in 
Section [TV] in Fig. [8] for U = 2 and U = 5 as a function 
of A. 

We find that the low value for n p h in the AFM state in- 
creases quite slowly as the coupling strength is increased. 
Once the transition to the CO state has occurred n p h 
increases substantially. The expectation values connect 
continuously, but increase rapidly for larger A. The value 
for the normal state lies below the CO result but increases 
rapidly for A ~ 3, where the metal-bipolaron transition 
occurs. This behavior can be compared with the result 
for larger U = 5, which is shown in Fig. [5] (bottom 
right). We can see that n p h, being small in the AFM 
state, remains nearly unaltered when the coupling is in- 
creased. The large U and strong AFM order suppresses 
local charge fluctuations. The value increases rapidly 
near U c s — in the CO state. The expectation value 
changes therefore discontinuously at the transition. The 
behavior in the normal state is similar to the earlier case 
and the metal-bipolaron transition (A ~ 5.3), which is 
also discontinuous for these parameters. 

The charge order, which we have characterized by $ co , 
can also be seen directly in the displacement expectation 
value on the A-sublattice (xa) — (bA + b' A ) / \J2uiq = (x) . 
This value is always zero in the AFM and N state but 
finite once the CO symmetry is broken. In Fig. [15] it is 
plotted for various values of U as a function of A. 
We can see that similar to the behavior of the order pa- 
rameter (x) increases close to the transition, and more 
rapidly for larger U. In its dependence on the coupling 
strength, it appears to be very similar to the order pa- 
rameter (Figs. [2] and 3]). In fact, one can show that they 
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FIG. 15: (Color online) The expectation values (a;) for various 
values of U as a function of A for ljq — 0.6 in the CO state. 



are directly related by an exact identity, 

2-y/A. 



(x) 



(9) 



which can be derived by considering an additional term 
H l c = u>oc(bi + b\) in the Hamiltonian and calculating the 
derivative with respect to c. The numerical values for the 
left and right hand side of equation (j9]) agree very well. 
We can see that the slope at the transition increases with 
U. At large U, or very small uiq, there is a very sharp 
transition as seen for the case U = 5. 

The effect of the strong electron lattice coupling can 
be seen in the displacement fluctuations (i 2 ), which are 
plotted in Fig. [16] for U — 2 and U = 5 and a range of A. 
The behavior is reminiscent of the phonon expectation 
value n p h, where a continuous rise near the transition is 
visible for weaker coupling and a discontinuity from the 
AFM to the CO state of the stronger coupling case. Here 
also the normal BP state possesses a larger value than 
the CO state. The close comparison with the results for 
n p h shows that for large coupling in the BP and CO state 
a;Q(l+4n p h)/2 gives a good fit. This result can be derived 
by taking = m t \, = 1, and performing a displaced 
oscillator transformation to new phonon operators, a, a\ 
(t) = b^+g/ujQ. The ground state |gs) then corresponds 
to the state o|gs) = 0, and in this state (b^b) = <? 2 /w 2 = 
n ph , and (x 2 ) = Xq(1 + 4g 2 /w 2 )/2, giving the required 
result. The state |gs) in the original basis corresponds 
to the coherent state J2 n "^rl n ) (^l n ) = n \ n ))i with 
a = ~ ^/ n ph, and the result can alternatively be derived 
by taking the expectation value of x 2 in this state. The 
decoupled oscillator state is an eigenstate \ n), but when 
strongly coupled to the electronic system the nature of 
the state changes to the coherent ground state due to the 
displacement of the oscillator. 

If we multiply (x 2 ) by Wq/2 we obtain the potential 
energy of the harmonic oscillator. The comparison with 
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FIG. 16: (Color online) The displacement fluctuations (x 2 ) 
for U = 2, 5 as a function of A for ljq = 0.6 for the AFM, CO 
and normal state. 



-Eph shows then that most of the phonon energy is in the 
potential energy due to the charge redistribution, and 
only small proportion in the kinetic energy of the oscil- 
lator. 

The real lattice fluctuations are large in the BP state 
where the local occupancy changes from double to zero, 
but (x) = 1. A measure of these fluctuations is the 
quantity Ax 2 = ((x - (x)) 2 } = (x 2 ) - (x) 2 . This is 
a much smaller quantity than (x 2 ) in the CO state and 
only large near the transition, as can be seen in Fig. [T7] 
In the uncoupled state Ax 2 /xq = (x 2 )/xq = 1/2. In the 
displaced oscillator (coherent) state, which describes the 
strong coupling situation, (x) 2 = 2x^g 2 /uJq, so combining 
this with the expression derived earlier for (x 2 ) , one again 
finds Ax 2 /xq = 1/2. Both limiting cases can be found in 
Fig. [T71 where Ax 2 increases with A in the AFM state. 
It then falls again to 1/2, when the system is strongly 
ordered. So the "lattice fluctuations" are largest at the 
transition. 

It is possible to use the density matrix approach in 
the NRG 50 together with the real space harmonic oscil- 
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FIG. 17: (Color online) The expectation values Ax 2 for var- 
ious values of U as a function of A for loq = 0.6 in the AFM 
state for A < U and in the CO state for A > U. 



lator eigenfunctions to compute the oscillator displace- 
ment probability function P(x), which gives further in- 
sight into the behavior of the phonons. Details for this 
method have been presented elsewhere^. The moments 
of this distribution function, (x m ) — Jdx P(x)x m , can be 

calculated from P(x) and are in agreement with the value 
determined from the groundstate expectation values. 



VII. SPECTRAL PROPERTIES 



In this section we turn to the excitation spectra of the 
coupled electron-phonon system. We first consider the 
electronic local lattice Green's function G a ^(Lu), which 
is given by the momentum sum of the diagonal element of 
(|2|), and its spectral function p a ^ a {uj) = — ImG a ,a{i>j)/Tr. 
In these calculations the sublattice self-energies Y, ata (oj) 
enter, which are calculated according to Eq. (|B1I) fol- 
lowing Bulla et ali£2. We will comment on complications 
which can arise in appendix [Bj Pa.ai^) is of special in- 
terest in studying the broken symmetry behavior. In the 
normal state the Green's functions are the same for the 
sublattices as well as the spin projections. For the cases 
with symmetry breaking, CO and AFM, these function 
will differ, i.e. for the CO case the sublattice Green's 
functions differ and for the AFM case the different spin 
projections. We focus on the A-sublattice majority spin 
spectral function pa,-\{u)- Note that at half filling the 
spectra for minority spin in the AFM case, pa,i(u), and 
for the B-lattice for the charge order, pb^(oj), can be 
obtained from ui — > — oj. In order to calculate the full 
electron Green's function one has to put the different 
sublattice Green's functions together—. 



A. Electron spectra 

We first consider the electronic spectral functions near 
the transition for the cases of U = 2 and U = 5. In Fig. 
1181 we plot the N state spectral function in comparison 
with the corresponding symmetry broken one. We have 
included a N state spectrum for A — > in order to first 
see the effect of the phonons and associated modification 
in the normal state. 
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FIG. 18: (Color online) The local A-lattice spectral functions 
in comparison for (7 = 2 (upper panel) and (7 = 5 (lower 
panel). Left: Comparison of N state spectrum with AFM 
majority spin pA,f(uj) near the transition. Right: Compar- 
ison of N state spectrum with CO A-site pA,t( tJ ) near the 
transition. 



For U = 2 in the upper part of the figure we see that close 
to the transition the N state spectrum deviates little from 
the A — > situation, with only a small extra renormal- 
ization of the low energy excitations [z(X = 0) = 0.73 
and z(X = 2) = 0.57]; the metal bipolaron transition 
with z — > occurs for a larger A ~ 3. The imaginary 
part of the electronic self-energy due to electron phonon 
scattering becomes finite when > ujq, where ujq is the 
renormalized phonon frequency (see Sec. IVH Bp . 

For the AFM and CO symmetry broken state we can 
see characteristics of the weak coupling instability at the 
Fermi surface {uj = 0) in the sublattice spectral function. 
It fits well to the mean field description, where a square 
root divergence is found below the gap^ The higher en- 
ergy parts are little modified for the case of U = 2 apart 
from the broadening of the band edges, but no features 
which can be attributed to the phonons can be identified. 

In the lower panel of Fig. [T5] we can see the situa- 
tion for U = 5. For A — > we have the well-known three 
peak structure with lower and upper Hubbard band and a 
quasiparticle peak in the N state. For A ~ 5 the quasipar- 
ticle structure is still visible, but the Hubbard peaks have 
been modified to high energy shoulders and cannot really 
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be recognized any more. The effect is as if the effective 
electron-electron interaction is screened by the phonons. 
The quasiparticle weight becomes larger in this regime 
on increasing A [z(A = 0) = 0.1 and z(X = 5) = 0.27, 
see also Fig. QT]. The low energy features of the spectra 
with symmetry breaking look similar. Directly above and 
below the spectral gap one sees pronounced peaks with 
larger weight for the one below the gap. At higher en- 
ergies the spectra look different in both situations. Here 
both features from the large U as well as from higher or- 
der polaronic behavior can play a role. When interpreting 
the spectra one has to take into account the broadening 
and the limited energy resolution of the NRG at higher 
energies, which limit the accuracy. The AFM state with 
strong electron phonon coupling seems to show polaronic 
behavior at multiples of wo with decreasing weight. The 
charge ordered spectrum shows a principal peak at a po- 
sition, which is a bit less than the fully polarized mean 
field shift Un^ a — 2A$ co w 2.5. The large couplings for 
this case play a role on different energy scales. One may 
note that for the given parameters the order parameters 
are rather large, $ w 0.4 (see Fig. whereas the spec- 
tral gap is relatively small. 

To study the behavior of the spectral functions at dif- 
ferent couplings we give various plots in Fig. 1191 On 
the left hand side we keep U = 2 fixed and show (from 
top to bottom) the spectra for N, AFM and CO state for 
various A. On the right hand side U = 5 is kept fixed. 

The top panel shows the metal to bipolaron transition 
when A is increased. One can see the strong narrowing 
of the quasiparticle band near the transition, which is 
accompanied by z — > 0. A spectral gap develops when 
A exceeds a critical coupling. The details of the spec- 
tral functions have been analyzed by Koller et alJ^. In 
the AFM case in the middle we can see for the weaker 
coupling case how, on increasing A, the AFM order and 
magnitude of the spectral gap decreases. The electron 
phonon coupling is effective here in screening the repul- 
sive [/-term. No polaronic features can be identified in 
the spectra as the coupling is fairly weak. At stronger 
coupling, the AFM state is hardly affected for a range of 
A. When approaching the transition we find visible mod- 
ifications of the spectral functions including a reduction 
of the spectral gap and polaronic peaks. 

In the bottom part the spectra in the CO state can be 
seen. Near the transition the spectra differ for the weak 
and strong coupling case as discussed before, but when 
A exceeds U by a certain amount the spectra look very 
similar, and the different U term is not directly visible 
anymore. As noted before the main peak for the sub- 
lattice spectra is located near the mean field shift and 
its position moves linearly with A as expected. There 
is a pronounced quasiparticle peak near the transition 
which becomes suppressed for larger values of A. This 
suppression can be partly due to the broadening in the 
NRG procedure as discussed in detail for superconduct- 
ing solutions^. 
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FIG. 19: (Color online) The local A-lattice spectral functions 
in comparison for (7 = 2 (left) and (7 = 5 (right) for various 
A and ujo — 0.6. Top panel: N state spectrum. Middle panel: 
AFM majority spin pA, t(w). Lower panel: CO A-site pA,t{w). 



B. Phonon spectra 

In this section we study the spectral properties of the 
phonons. From these we can find out how the excitations 
of the bosonic sector are modified through the interaction 
with the electronic system. Especially near a transition 
a strong phonon softening can be indicative of a lattice 
modification or instability. We consider the function 



(10) 



which can be calculated in the NRG from the matrix 
elements and excitations. The spectral function, Pb{u) = 
—lmB(u))/Tr, has the properties at T = 0, 



dw pb(<jj) = 1, /da; /?t,(w) 



''phi 



and the free propagator has the form 

B°M - -r^ ■ 



— UJq 



(11) 



(12) 



We define a phonon self-energy £ p h(w) the full propaga- 
tor reads 

BH = — ~ 1 ^ — ■ (13) 



- £ ph (u;) 
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For a decoupled electron-phonon system S p h(w) = 
such that p&(w) is a delta function peaked at w . In 
the strongly interacting system the mode can be renor- 
malized to and broadened, however, no q- dependence 
develops in the infinite dimensional model. 

Again we focus on the cases with fixed U = 2 and 
U = 5 for variable A, and we compare the results from 
the N state with the ordered state. In Fig. [30] we plot 
in the upper part the results for the spectral function of 
the phonons pb{ui) for U — 2 for the N state (left) and 
the ordered states (right) and the lower part for U = 5. 

When the electron and phonon systems are weakly 
coupled the expected delta- function for Pb{u) is found. 
When A is increased the phonon mode is substantially 
renormalized and softens, most markedly at the M-BP 
transition for A ~ 3 where ujq —¥ 0. The negative spec- 
tral weight for co, which builds up there, is directly related 
to the phonon expectation value rt p h- In the BP state it 
is not resolved any more. The phonon mode then hard- 
ens back to uq. Apart from the softening we also find a 
broadening of the phonon spectrum, when the system is 
strongly coupled. The behavior for U = 5 is similar to 
the U = 2. A more detailed discussion can be found in 
Ref. HE 
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FIG. 20: (Color online) The local phonon spectral functions 
Pb{<jj) in comparison for U — 2 (upper panel) and U = 5 (lower 
panel). Left: N state. Right: AFM state (smaller values of A 
and CO state for the larger values of A. 

In the right hand side panels of Fig. [201 the correspond- 
ing behavior of /9&(w) is shown in the states with broken 
symmetry. The results for A < U correspond to the AFM 
phase and the ones for A > U to the CO phase. Com- 
paring the metal-bipolaron transition to the AFM-CO 
transition we find that the effect of the phonon softening 
is much reduced. For U — 2 we can see a visible effect 
that the oscillator mode is renormalized to ojq ~ 0.27 
and broadened. The effect is comparable to the normal 
state for the same values of U and A, but not as strong 
as at the M-BP transition. The softening is of similar 



magnitude in the stronger coupling case [7 = 5. There 
a significant broadening of the phonon mode visible near 
the transition due to the large coupling of the phonons 
to the electronic system. When the system is well in the 
ordered state, AFM or CO, the phonon dynamics is lit- 
tle modified by the electronic system and a nearly free 
phonon mode is observed both for large U 3> A in the 
AFM state and U <C A in the CO state. The sum rules 
for Pb(oj) given in equation (jlip are satisfied to within a 
few percent for smaller values of A. For larger A in the 
BP and CO state the sum rules are not well satisfied due 
to reasons discussed in Ref. Uji 



VIII. CONCLUSIONS 

In this study of the competing interactions in the 
Holstein-Hubbard model we have examined the transi- 
tions to both AFM and CO, both in the weak and strong 
coupling regimes. To lowest order, the effective frequency 
dependent interaction U e s{cj) between the electrons is 
given by 



C/effM =U + 



2g 2 u 



(14) 



the second retarded term arising from phonon exchange. 
On the lowest energy scale to = 0, U e g(0) — U — A 
(A = 2g 2 /ujo), and the sign of this interaction depends 
on the relative strength of U and A. In studying the 
competition between AFM and CO, it is not surprising 
to find that the transition between these states occurs 
when U - A, as U eS (0) > favors AFM and U eS (0) < 0, 
the CO state. What is a surprising result of this study 
is that this condition still has some validity in the strong 
coupling regime, when both U and A are large, as the 
transition is still found to occur when U ~ A. We find, 
however, the nature of the transition does depend on 
the strength of the couplings, and also the phonon fre- 
quency ujq. The transition is found to be continuous for 
weak couplings, and a high phonon frequency ujq, but be- 
comes discontinuous in the strong coupling regime, and 
for smaller values of coq. 

To gain further insight into this result, we have looked 
in detail at the quasiparticle excitations in the normal 
state. We have calculated both the quasiparticle weights 
z and the effective local quasiparticle interaction U r . We 
find that the local quasiparticle interaction U r changes 
sign when U ~ A, just in the region where the AFM-CO 
transition occurs; this is consistent with the interpreta- 
tion of the transition as due to a Fermi liquid instability. 
Though the interaction between the quasiparticle goes to 
zero at U ~ A, the quasiparticles may still be quite sig- 
nificantly renormalized. For example, for U = A = 5 we 
find z ~ 0.3. The fact that the local quasiparticle inter- 
action goes to zero in the region U ~ A, corresponding 
to ?7eff(0) = 0, suggests that the two terms contribut- 
ing to U e g(ui), are renormalized on the very low energy 
scale in a similar way. This is somewhat surprising, as in 
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considering similar competing interactions in the case of 
superconductivity, it is generally assumed that the dom- 
inant renormalization is of the Coulomb term so that it 
does not overwhelm the attractive term from phonon ex- 
change. That the two terms are renormalized here in a 
similar way may be a feature of the Holstein-Hubbard 
model, where the phonon term is coupled to the occupa- 
tion of a local charge. A model in which the phonons are 
coupled to a redistribution of the local charge, as with a 
coupling to Jahn- Teller modes, might behave differently. 
This topic deserves further investigation. 

In calculating the individual contributions to the to- 
tal energy in the different ordered states and the normal 
state, we have been able to show the subtle interplay of 
the various terms. These vary in the weak and strong 
coupling regimes, and may change discontinuously at the 
transition. They also depend on the phonon frequency 
ojq. In the weak coupling regime the energy gain in the 
broken symmetry state is via a reduction of the potential 
energy relative to that of the normal state, whereas at 
strong coupling it is the kinetic energy which is lower in 
the ordered state. This appears to be a general feature. 

The final part of this study has been concerned with 
the spectra, both of the electrons and phonons. The main 
effects seen in the phonon spectra are a softening and a 
broadening of the phonon mode in the region of the tran- 
sition. In the AFM and CO states well away from the 
transition there is little effect of the coupling to the elec- 
trons on the phonon spectrum. It is more difficult to 
summarize the results for the electron spectra, as there 
significant differences develop on all energy scales as the 
interaction parameters are varied, and as the long range 
AFM or CO develops. For relatively small values of U 
and A the main differences are in the region near the 
Fermi level for AFM or CO states compared with the 
normal state. This is due to the development of the sub- 
lattice structure. For large values of U and small values of 
A there is the triple peak structure of the Hubbard model, 
with the narrow renormalized quasiparticle band at the 
Fermi level flanked by the broadened 'atomic' peaks. As 
A is increased to A ~ U , in the normal state, the narrow 
quasiparticle band persists, though broadened somewhat, 
and the atomic-like peaks broadened into shoulders. In 
the AFM or CO states, the quasiparticle band at the 
Fermi level develops the features associated with the sub- 
lattice structure, as in the weak coupling case. This is 
also accompanied by much larger shifts of spectral weight 
on the high energy scales in the sublattice spectral den- 
sity. 
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Appendix A: Mean field theory in the adiabatic limit 

For the mean field theory in the adiabatic limit, the 
starting point is the Hamiltonian in the form 



(Al) 



2M + 2 



+5F Y %i ( H ^i," ~ *) + H 
i (7 i 

where the parameters of (|A1[) and (UJ are related by 
ujq = y/k/M, <7f = V^og, and A = g%/k. In this 
Hamiltonian we can take the limit M — > oo, such that 
the kinetic term for the phonons vanishes and we replace 
the operator ii by a static field a;,-. From this we obtain 
in mean field theory the potential 

k 

V(xi) = Y ~ x i + 9f Y Xi ( n * - X ) + ^ kin + E v> ( A2 ) 



where = ^ 
dV(xi)/dxi — yields, 



The condition for a local minimum 



— (rii-1). 



(A3) 



We restrict ourselves to homogeneous solutions, and from 
the Hamiltonian the mean field self-energy can be read 
off, 



£ Q , CT H = tM" CT -A(n Q -l), 



(A4) 



independent of u and we have employed (|A3[) for the 
second term. The index a — A, B corresponds to the 
sublattice and a to the spin. 

In order to determine = (n™) we need to consider 
the equation 

oo 

H " = JjJ2( C Lk,a C »,k,a) = JjY dUJ /-M/Vfc.aH, 



where 



C4,<t(w + )Cb,ct(w H 



(A5) 



with uj + = lj + if}. We have used the matrix Green's 
function for the bipartite lattice in the form ©. This is 
most easily evaluated with the identity 



Tm- 



(a,a(u + )/K 



CaAu + )(bA" + )-4 



W m,<r( £ fe) (5 [ W - W m,a(£fc)]- 



m— ± 



The excitation u) m ,<r(£k) are determined from the poles 
of the Green's function, 



CA, CT (w)CB,a(w) -4 = 0, 



(A7) 
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which yields generally 
w±,o-(£fe) = ± ^fc, CT , (A8) 



where 



From this we can see that if the order parameters are 
equal and exclusive the CO state has lower energy for 
A > U and the AFM state otherwise. 

Appendix B: Calculation of the self-energy 



&k,o = ye h + — ■ . 



(A9) 



The weights a (ek) are generally given by the inverse 
of the derivative w.r.t oj of 



In NRG calculations it is common practice to deter- 
mine the self-energy from the Green's function G a _ a (uj) 
and the higher order Green's function F Qj(T (u;) via^ 



(A10) 



evaluated at uj± t a-{£k), 



u m,A £k > ~ 2 i ? 7 7 \\2- { AU ) 

Using these result we find 



(A12) 



through which can be determined self-consistently. 

Once is determined we can calculate the ground 
state energy to determine which state has the lowest en- 
ergy. The expression for the total energy reads 

Emi = M^2( £h ( C A,k,a C B,k,a) + h.C.) 



k.a 



where we have substituted (|A3[) for a;. This can also be 
written as 

E mi = ^-A[(^ ) 2 + ($f ) 2 ] 

where 



pmt _ \ - f, Po(e)g 2 ( 1 



For half filling, $ co = = $ afm = |**J 



I "^afm I this can be wr itten in the simple form 



EU = E^ n + (U- 2A)$ 2 D - C/$ 2 fm + | (A13) 



S q ,<t(w) = C/ 



(Bl) 



This can be derived in an equations of motion approach. 
As F, G are complex functions, F = F R + iF l , we can 
write 



U 



f r g r + F i G i + ,( F i G R _ F R g I) 
(G R ) 2 + (G 1 ) 2 ' 



(B2) 



where we have omitted the indices and the arguments. 
The procedure (|B1|) for obtaining £ has turned out to 
work well in many cases both for impurity models and 
lattice models within the DMFT framework, ?§ The imag- 
inary part of the retarded self-energy has the well-known 
property ImE a (T (w) < 0, which is respected in Eq. (|B2[) 
if F and G are the exact Green's functions. However, 
in a numerical self-consistent DMFT calculation of G, F 
small inaccuracies - usually near = - can lead to 
F i (uj)G r (u)-F k (u))G 1 (u) > and thus slightly positive 
values for ImE via Eq. (|B2[) . Clearly this is physically 
incorrect. We have used two different ad-hoc procedures 
to deal with this complication. The first one (a) is to sub- 
tract from a (oS) the values by which it exceeds zero 
in a certain interval around uj = 0. The second one (b) 
is to cut-off Ejj, a {ijo) at zero, i.e. to set it equal to zero 
for all values where it is positive. We found that in most 
cases the procedures give approximately the same result. 
However, very close to the AFM-CO transition the pro- 
cedure can have an effect on the final result obtained via 
the self-consistency equation. One finds that ordered so- 
lutions are a bit less stable for method (a). We have 
decided to present all results in this paper obtained by 
using method (b). 
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